4. Iterative Solvers

So far we have used direct solvers to solve the linear system of equations. Although a direct solver can profit from the sparse matrix, it’s arithmetic complexity is sub-optimal (i.e. not linear in the number of degrees of freedom). For large-scale problems iterative solvers are a must.

The conjugate gradient (cg) method is the standard method for symmetric and positive definite matrices. It’s convergence rate depends on a preconditioner, what is a cheap approximative inverse to the matrix.

from ngsolve import *
from ngsolve.webgui import Draw

We generate a 3D geometry and mesh using the OCC constructive solid geometry (CSG) modeler:

from netgen.occ import *
cube = Box((0,0,0),(1,1,1))
cyl = Cylinder((0,0.5,0.5),X, r=0.2, h=1)
cube.faces.name = "outer"
cyl.faces.name = "cyl"
shape = cube-cyl
Draw(shape);

Generate a mesh, and perform uniform mesh refinement:

ngmesh = OCCGeometry(shape).GenerateMesh(maxh=0.1)
for l in range(1):
    ngmesh.Refine()
mesh = Mesh(ngmesh)
mesh.Curve(3)
Draw (mesh);

fes = H1(mesh, order=3, dirichlet="outer", wb_withedges=False)
print ("we have", fes.ndof, "unknowns")
u,v = fes.TnT()

a = BilinearForm(grad(u)*grad(v)*dx)
f = LinearForm(v*dx)

# c = Preconditioner(a, "direct", inverse="sparsecholesky")
c = Preconditioner(a, "local")
# c = Preconditioner(a, "bddc")

gfu = GridFunction(fes)

we have 177430 unknowns

assemble system and setup preconditioner in parallel:

with TaskManager():
    a.Assemble()
    f.Assemble()

solve the system using the preconditioned conjugate gradient method:

from ngsolve.krylovspace import CGSolver

with TaskManager():
    inv = CGSolver(mat=a.mat, pre=c.mat, printrates=True, maxiter=400)
    gfu.vec.data = inv * f.vec

CG iteration 1, residual = 0.024273223974312072     
CG iteration 2, residual = 0.03720447589534199     
CG iteration 3, residual = 0.04829608173279241     
CG iteration 4, residual = 0.04783145202720858     
CG iteration 5, residual = 0.03503614777663428     
CG iteration 6, residual = 0.02937808518603648     
CG iteration 7, residual = 0.031006371846374843     
CG iteration 8, residual = 0.02881189922538989     
CG iteration 9, residual = 0.026991393488043684     
CG iteration 10, residual = 0.02560781323347911     
CG iteration 11, residual = 0.023059227507159427     
CG iteration 12, residual = 0.0201887142527063     
CG iteration 13, residual = 0.017677856690926363     
CG iteration 14, residual = 0.01570390414750895     
CG iteration 15, residual = 0.014508920088754975     
CG iteration 16, residual = 0.013060236654630368     
CG iteration 17, residual = 0.011487172720698016     
CG iteration 18, residual = 0.010359419640115412     

CG iteration 19, residual = 0.009639280718573355     
CG iteration 20, residual = 0.008983784519869012     
CG iteration 21, residual = 0.00825049617120831     
CG iteration 22, residual = 0.007352872683189113     
CG iteration 23, residual = 0.006439789142169032     
CG iteration 24, residual = 0.005596514201958335     
CG iteration 25, residual = 0.0047662198736664165     
CG iteration 26, residual = 0.004022073286663706     
CG iteration 27, residual = 0.0034196513449910827     
CG iteration 28, residual = 0.0029400048222059856     
CG iteration 29, residual = 0.0025083332705242     
CG iteration 30, residual = 0.002104649927409853     
CG iteration 31, residual = 0.0017484379251981059     
CG iteration 32, residual = 0.0014397425850613986     
CG iteration 33, residual = 0.0011912144262735038     
CG iteration 34, residual = 0.0009893154642094015     
CG iteration 35, residual = 0.0008329679456139164     
CG iteration 36, residual = 0.000727286112694378     
CG iteration 37, residual = 0.0006352763333349189     
CG iteration 38, residual = 0.0005645223644237872     
CG iteration 39, residual = 0.0004942564003676568     
CG iteration 40, residual = 0.00042780257753367347     
CG iteration 41, residual = 0.00036683687047396385     
CG iteration 42, residual = 0.0003112297948991815     
CG iteration 43, residual = 0.00026555485540291635     
CG iteration 44, residual = 0.00023057965122062087     
CG iteration 45, residual = 0.00020810720521774679     
CG iteration 46, residual = 0.000192089094018006     
CG iteration 47, residual = 0.0001785379289683777     
CG iteration 48, residual = 0.00016491550344460956     
CG iteration 49, residual = 0.00014983626453891514     
CG iteration 50, residual = 0.00013386486143320445     
CG iteration 51, residual = 0.00011598796595447697     
CG iteration 52, residual = 9.959344694700464e-05     
CG iteration 53, residual = 8.519641137784835e-05     
CG iteration 54, residual = 7.322459830883028e-05     
CG iteration 55, residual = 6.417417797477794e-05     
CG iteration 56, residual = 5.69968450773374e-05     
CG iteration 57, residual = 5.166607912750661e-05     
CG iteration 58, residual = 4.711071719270871e-05     
CG iteration 59, residual = 4.323161539942034e-05     
CG iteration 60, residual = 3.952383801108806e-05     
CG iteration 61, residual = 3.548362305135951e-05     
CG iteration 62, residual = 3.115776481808801e-05     
CG iteration 63, residual = 2.6858896178966227e-05     
CG iteration 64, residual = 2.2452143271619636e-05     
CG iteration 65, residual = 1.8667683822022986e-05     
CG iteration 66, residual = 1.5676517902103544e-05     
CG iteration 67, residual = 1.3491991701423917e-05     
CG iteration 68, residual = 1.1932267987396299e-05     
CG iteration 69, residual = 1.1006946297625633e-05     
CG iteration 70, residual = 1.0456629339641562e-05     
CG iteration 71, residual = 9.883740543600476e-06     
CG iteration 72, residual = 9.186544775446457e-06     
CG iteration 73, residual = 8.152829147699973e-06     
CG iteration 74, residual = 6.959623023961479e-06     
CG iteration 75, residual = 5.777798497628985e-06     
CG iteration 76, residual = 4.77500717595072e-06     
CG iteration 77, residual = 3.979591970804925e-06     
CG iteration 78, residual = 3.3697656756024886e-06     
CG iteration 79, residual = 2.9314682049183286e-06     
CG iteration 80, residual = 2.6063431812845677e-06     
CG iteration 81, residual = 2.3391092581870584e-06     
CG iteration 82, residual = 2.107326811990917e-06     
CG iteration 83, residual = 1.878388180936518e-06     
CG iteration 84, residual = 1.6729158041834687e-06     
CG iteration 85, residual = 1.4773593883411322e-06     
CG iteration 86, residual = 1.2787374798880876e-06     
CG iteration 87, residual = 1.1065178720984112e-06     
CG iteration 88, residual = 9.440248921041147e-07     
CG iteration 89, residual = 8.07063285821953e-07     
CG iteration 90, residual = 6.801607317244886e-07     
CG iteration 91, residual = 5.689337874208956e-07     
CG iteration 92, residual = 4.77978366369197e-07     
CG iteration 93, residual = 4.032311010598506e-07     
CG iteration 94, residual = 3.443374548684612e-07     
CG iteration 95, residual = 2.984216523076836e-07     
CG iteration 96, residual = 2.65037355337294e-07     
CG iteration 97, residual = 2.3871960758477033e-07     

CG iteration 98, residual = 2.1411228208313604e-07     
CG iteration 99, residual = 1.9221188353254465e-07     
CG iteration 100, residual = 1.7291785140323965e-07     
CG iteration 101, residual = 1.538520444575197e-07     
CG iteration 102, residual = 1.343657564899277e-07     
CG iteration 103, residual = 1.1423742371249511e-07     
CG iteration 104, residual = 9.671589404334885e-08     
CG iteration 105, residual = 8.178565207017999e-08     
CG iteration 106, residual = 6.952013932928278e-08     
CG iteration 107, residual = 5.900821337943827e-08     
CG iteration 108, residual = 5.0648361308266256e-08     
CG iteration 109, residual = 4.39992418314956e-08     
CG iteration 110, residual = 3.8905133069483565e-08     
CG iteration 111, residual = 3.525703762831813e-08     
CG iteration 112, residual = 3.252791286211454e-08     
CG iteration 113, residual = 2.966916187816114e-08     
CG iteration 114, residual = 2.6685243081537716e-08     
CG iteration 115, residual = 2.3725137783681394e-08     

CG iteration 116, residual = 2.113450100005827e-08     
CG iteration 117, residual = 1.849349701822083e-08     
CG iteration 118, residual = 1.6151023745212e-08     
CG iteration 119, residual = 1.3930393942545725e-08     
CG iteration 120, residual = 1.1937000391384313e-08     
CG iteration 121, residual = 1.0319921823881497e-08     
CG iteration 122, residual = 9.056959370552789e-09     
CG iteration 123, residual = 8.008991901523645e-09     
CG iteration 124, residual = 7.154633406243156e-09     
CG iteration 125, residual = 6.394790148219618e-09     
CG iteration 126, residual = 5.746953924346185e-09     
CG iteration 127, residual = 5.180088590345217e-09     
CG iteration 128, residual = 4.62984083620587e-09     
CG iteration 129, residual = 4.15123841287872e-09     
CG iteration 130, residual = 3.717167939488051e-09     
CG iteration 131, residual = 3.2863666418280675e-09     
CG iteration 132, residual = 2.875922737363393e-09     
CG iteration 133, residual = 2.4976440982222984e-09     
CG iteration 134, residual = 2.169561335323967e-09     
CG iteration 135, residual = 1.868479538298537e-09     
CG iteration 136, residual = 1.618615763529619e-09     
CG iteration 137, residual = 1.396560493915583e-09     
CG iteration 138, residual = 1.2300628068545116e-09     
CG iteration 139, residual = 1.083725370099002e-09     
CG iteration 140, residual = 9.72876629737617e-10     
CG iteration 141, residual = 8.73270816130791e-10     
CG iteration 142, residual = 7.863948798897759e-10     
CG iteration 143, residual = 6.947907452242882e-10     
CG iteration 144, residual = 6.079255156456969e-10     
CG iteration 145, residual = 5.246752600351578e-10     
CG iteration 146, residual = 4.5036587153538776e-10     
CG iteration 147, residual = 3.8833470765450594e-10     
CG iteration 148, residual = 3.3461037938448853e-10     
CG iteration 149, residual = 2.922388858222119e-10     
CG iteration 150, residual = 2.5493981049920296e-10     
CG iteration 151, residual = 2.2453025075986217e-10     
CG iteration 152, residual = 1.978682508573511e-10     
CG iteration 153, residual = 1.7309758865948704e-10     
CG iteration 154, residual = 1.5142752452582156e-10     
CG iteration 155, residual = 1.3191628519621144e-10     
CG iteration 156, residual = 1.1577503459937672e-10     
CG iteration 157, residual = 1.0234668879908436e-10     
CG iteration 158, residual = 9.025208464709837e-11     
CG iteration 159, residual = 8.016629210425645e-11     
CG iteration 160, residual = 7.059681741993443e-11     
CG iteration 161, residual = 6.18138419616014e-11     
CG iteration 162, residual = 5.4119013093614436e-11     
CG iteration 163, residual = 4.6976988867082243e-11     
CG iteration 164, residual = 4.060099093133944e-11     
CG iteration 165, residual = 3.483496350815462e-11     
CG iteration 166, residual = 2.96277528769076e-11     
CG iteration 167, residual = 2.511684953302571e-11     
CG iteration 168, residual = 2.1314124458049647e-11     
CG iteration 169, residual = 1.839012228744114e-11     
CG iteration 170, residual = 1.5920034082411392e-11     
CG iteration 171, residual = 1.3931214825016162e-11     
CG iteration 172, residual = 1.224541437575706e-11     
CG iteration 173, residual = 1.0842719389574334e-11     
CG iteration 174, residual = 9.609156129884664e-12     
CG iteration 175, residual = 8.567563010304951e-12     
CG iteration 176, residual = 7.63605044217954e-12     
CG iteration 177, residual = 6.77241600237929e-12     
CG iteration 178, residual = 5.745193864781057e-12     
CG iteration 179, residual = 4.846940156879799e-12     
CG iteration 180, residual = 4.121887000467237e-12     
CG iteration 181, residual = 3.4924718882163277e-12     
CG iteration 182, residual = 2.946742557391849e-12     
CG iteration 183, residual = 2.4886293668436434e-12     
CG iteration 184, residual = 2.1226341228537092e-12     
CG iteration 185, residual = 1.8325206620203842e-12     
CG iteration 186, residual = 1.5897837501020638e-12     
CG iteration 187, residual = 1.4001146028477738e-12     
CG iteration 188, residual = 1.2347868852528434e-12     
CG iteration 189, residual = 1.0935945240239527e-12     
CG iteration 190, residual = 9.777037400031674e-13     
CG iteration 191, residual = 8.798918129380025e-13     
CG iteration 192, residual = 7.825758317382021e-13     
CG iteration 193, residual = 6.789726755673787e-13     

CG iteration 194, residual = 5.798382194518679e-13     
CG iteration 195, residual = 4.859229251905611e-13     
CG iteration 196, residual = 4.051387668153864e-13     
CG iteration 197, residual = 3.390456762084308e-13     
CG iteration 198, residual = 2.828305254963196e-13     
CG iteration 199, residual = 2.3546767272707774e-13     
CG iteration 200, residual = 1.9849692018148548e-13     
CG iteration 201, residual = 1.7126408342545894e-13     
CG iteration 202, residual = 1.504697886729289e-13     
CG iteration 203, residual = 1.3443993743744704e-13     
CG iteration 204, residual = 1.218333246195873e-13     
CG iteration 205, residual = 1.1113244400782543e-13     
CG iteration 206, residual = 1.0098118188191275e-13     
CG iteration 207, residual = 8.978667420785926e-14     
CG iteration 208, residual = 7.789492711962835e-14     
CG iteration 209, residual = 6.608271736967291e-14     
CG iteration 210, residual = 5.5217390371954335e-14     
CG iteration 211, residual = 4.5618798522729295e-14     

CG iteration 212, residual = 3.783428307207398e-14     
CG iteration 213, residual = 3.165288846010291e-14     
CG iteration 214, residual = 2.6689455992858683e-14     
CG iteration 215, residual = 2.331336429160311e-14     

Draw (gfu, draw_vol=True);

Exercise:

• Experiment with differ problem sizes and preconditioners

• How big systems can you solve on your computer (watch out for memory usage)